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Ubi  OtSTRMUTION  COOf 


Measurements  of  droplet  dispersion  and  vaporization  rates  are  carried  out  and  compared 
with  computational  predictions.  Droplet  dispersion,  mean  axial  velocity  and  the  axial 
autocorrelation  function  of  the  droplets  are  measured  using  a  laser  sheet  technique.  The 
experimental  facility  is  modified  to  accommodate  a  true  spray  using  an  ultrasonic  atomizer 
and  droplets  injected  with  a  fluorescent  dye.  The  computational  part  of  the  project  is  the 
development  of  a  Large  Eddy  Simulation  method  for  the  turbulent  jet.  A  LES-model  for 
the  unresolved  scales  of  the  flow  is  implemented  into  the  Navier-Stokes  solver.  The  results 
for  filtered  and  unfiltered  variables  are  compared  to  the  experiments  and  improvement  of 
the  particle  statistics  due  to  the  LES-model  is  found. 
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Research  Objectives 

The  overall  objective  of  the  project  is  to  test  current  modelling  of  droplet  scale  processes 
for  spray  combustion  in  a  well  defined  turbulent  shear  flow  by  compari.son  of  measurements 
of  droplet  dispersion  and  vaporization  rates  with  computational  predictions.  Our  specific 
goals  in  this  year  of  the  project  were  to  measure  the  autocorrelation  function  of  the  axial 
velocity  component  of  the  droplets  and  to  compare  it  with  stochastic  simulation  re.sults 
and  the  modification  of  the  experimental  facility  to  accommodate  a  true  spray  using  an 
ultrasonic  atomizer  and  droplets  injected  with  a  fluorescent  dye. 

The  computational  phase  of  the  project  is  aimed  at  the  development  of  a  Large  Eddy 
Simulation  method  for  the  turbulent  jet  that  has  been  the  object  of  the  experimental 
study.  The  goal  is  to  obtain  a  simtilation  of  the  jet  under  conditions  that  match  the 
experimental  conditions  so  that  droplet  di.spersion  and  vaporization  rates  can  be  compared 
and  correlations  of  drag  and  mass  transfer  evaluated.  Our  specific  goal  in  this  year  of  the 
project  was  to  implement  a  LES-model  for  the  unresolved  scales  of  the  flow  into  the  Navier- 
Stokes  .solver  and  to  compare  the  filtered  and  unfiltered  results  to  eva.luate  the  effect  on 
particle  dispersion  and  vaporization  rates. 
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Research  Accomplishments 


Experimental 


One  of  the  major  unknowns  in  modeling  particle  dispersion  in  turbulent  flows  is  the  form 
of  the  Lagrangian  velocity  autocorrelation  function  for  fluid  and  discrete  particles.  If  this 
function  is  known,  then  it  is  jio.s.sihle  to  obtain  the  dispersion  i.e.,  the  mean  square  displace¬ 
ment  from  a  known  initial  point.  Snyder  and  Lumley  (1971)  obtained  the  autocorrelation 
function  in  a  grid-generated  turbulence  but  there  are  no  data  for  particles  in  a  turbulent 
shear  flow  such  as  a  jet.  One  phase  of  the  project  over  the  last  year  has  been  devoted  to 
obtaining  this  fundamental  information. 


The  basic  experimental  methodology  has  been  reported  elsewhere  (Call  and  Kennedy, 
1991,  1992).  A  laser  sheet  is  formed  with  cylindrical  lenses.  Droplets  are  issued  from  a 
piezoelectric  generator  into  the  jet  near  the  centerline  of  the  jet  pipe  exit.  As  a  droplet 
passes  through  the  laser  sheet,  the  scattered  radiation  is  collected  and  imaged  onto  a 
position-sensing  photomultiplier  tube.  This  scheme  has  been  modified  to  obtain  velocity 
and  correlation  information.  Two  sheets  of  parallel  light  are  produced  with  a  beam  splitter 
and  right  angle  prism.  The  parallel  sheets  are  reflected  back  across  the  jet  with  a  retro- 
reflector  with  an  adjustal)le  spacing  between  the  pairs  of  laser  sheets.  The  time  between 
scattering  spikes  from  each  sheet  yields  droplet  velocity  data,  whatever  the  droplet  location 
across  the  jet  cross-section.  Autocorrelations  of  the  axial  particle  velocity  can  be  obtained 
from  the  two  velocity  measurements  from  the  two  pairs  of  sheets.  The  autocorrelations 
can  be  reported  in  either  an  Eulerian  form  as  a  function  of  spatial  separation  or  in  a  quasi- 
Lagrangian  form  as  a  function  c)f  time  from  the  first  transit  through  the  first  set  of  laser 
sheets. 


We  have  studied  hexadecane  droplets  ranging  in  diameter  from  50  pm  to  35  pm.  In 
addition,  hollow  glass  s]>heres  of  40  //m  diameter  have  been  used  to  approximate  fluid 
particles.  Figure  1  shows  the  quasi-Lagrangian  autocorrelation  function  for  axial  particle 


velocity  at  a  number  of  axial  locations  along  the  jet.  The  correlations  approximate  expo-  _  ,4*. 


nential  functions  quite  well.  Insufficient  data  are  available  at  longer  times  to  comment  on 
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the  applicability  of  other  functions,  such  a-s  the  Frenkiel  function.  The  Stokes  number  of 
the  particles  heis  a  clear  impact  on  the  correlation  function.  The  turbulence  Stokes  number 
(based  on  a  fluid  integral  length  scale)  is  0.16  and  0.33  for  the  35  and  50  micron  droplets 
respectively. 

Comparisons  \vfith  a  stochastic  simulation  have  emphasized  the  difficulty  inherent  in 
the  prescription  of  length  and  time  scales  in  turbulent  flows  with  a  cne- point  closure  model. 
The  stochastic  simulation  under  predicts  the  integral  length  scale  for  droplet  motion  by  a 
factor  of  two.  Of  course,  the  model  constant  that  prescribes  the  eddy  life  time  could  be 
adjusted  but  there  is  little  a  priori  guidance  for  the  choice  of  this  value. 

Droplet  dispersion  nieasureinents  have  been  carried  out  dining  the  last  year  with 
varying  jet  Reynolds  numlxn.  The  Reynolds  number  is  changed  by  adjusting  the  nozzle 
diameter  and  the  jet  exit  velocity  to  study  the  effect  of  the  turbulent  length  and  time 
scales  on  droplet  dispersion.  Typical  Reynolds  numbers  are  Re  —  10000,  20000,  and 
30000  obtained  with  the  nozzle  diameters  D  —  7mm,  10mm,  and  12.6mm.  Figure  1 
shows  the  average  axial  \-eIocify  and  the  root  mean  square  value  of  the  axial  velocity 
fluctuations  at  the  three  different  Reynolds  numbers.  The  turbulent  flow  at  the  higher 
Reynolds  numbers  will  also  create  noticable  radial  velocities  for  the  larger  (90pm )  droplets 
making  its  measurement  possible  and  facilitating  the  study  of  the  effect  of  the  turbulent 
scales  on  their  dispersion. 

The  experimental  facility  is  being  modified  to  accommodate  a  true  spray.  Figure  2 
shows  the  experimental  set  tip  for  the  spray  simulation.  An  ultrasonic  atomizer  has  been 
installed  in  a  chamber  to  provide  a  fine  spray  with  droplet  diameters  from  20  to  80  pm. 
The  ultrasonic  atomizer  offers  the  advantage  of  very  low  air  flow  rates  so  that  the  flow  of 
air  in  the  experiment  is  not  greatly  disturbed.  Single  droplets  that  contain  a  fluorescent 
dye  will  be  injected  by  the  piezoelectric  generator  onto  the  centerline  of  the  spray.  A 
holographic  edge  filter  will  remove  the  Mie  scattering  from  the  spray.  The  fluorescence 
from  the  dye- containing  droplet  will  be  detected  by  the  position-sensing  photomultiplier 
tube  as  before.  We  plan  to  study  the  impact  of  dispersed  phase  loading  in  the  spray  on 
droplet  dispersion  with  this  system. 
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Theory  and  Comvutation. 


Direct  simulation  of  turbulent  round  jets  is  still  not  within  the  capabilities  of  present  com¬ 
puters  (Reynolds,  1990).  Hence,  Large  Eddy  Sunulation  of  the  turbulence  in  round  jets  is 
the  most  realistic  approach  for  the  prediction  of  the  flow  field  and  truly  Lagrangian  particle 
dynamics.  The  simulation  of  turbulent  flow  fields  in  round  jets  is  based  on  accurate  finite- 
diffeience  methods,  which  offer  the  flexibility  necessary  for  the  treatment  of  non-periodic 
jet  flows  emitting  from  nozzles  and  the  consideration  of  a  variety  of  exit  conditions.  This 
aspect  of  the  project  was  desciiljed  in  detail  in  the  final  report  for  grant  AFOSR  89-0392 
and  further  details  can  be  found  in  Lienau.  Kennedy  and  Kollmann  ( 1993)  and  Lienau  and 
Kollmann  (1993)  of  the  ptiblication  list.  Two  new  contributions  were  accomplished  during 
the  period  1992-3.  The  numerical  treatment  of  the  coordinate  axis  r  =  0  for  unsteady  flows 
without  symmetries  was  analyzed  (details  in  the  appendix)  and  a  satisfactory  method  was 
found  to  avoid  the  loss  of  accuracy  near  the  axis.  The  second  and  main  contribution  was 
the  implementation  of  a  LES-model  for  the  non-resolved  scales  of  the  turbulence.  The 
Navier-Stokes  equations  are  written  for  filtered  variables 


fU.t) 


j  dx!G{x_  —  £  J) 


(1) 


where  (j(t  —  t',  f )  denotes  the  filter  function  and  /(x,  t)  a  dependent  variable.  The  Navier- 
Stokes  equations  (in  Cartesian  coordinates  for  convenience)  are  filtered  and  emerge  as 


di'c 

dxr 


=  0 
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The  modified  filtered  presstire  is  defined  by 


P=p+  ■^P^aPV'al’'3 


(4) 


and  the  correlations  of  filtered  and  sub-grid-scale  motion  ( v'^  =  Uq  —  Va )  are  given  by 


^aP  — 


(3) 
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(6) 


CoJ  =  I'al'J  +  VaV'., 


RnJ  S  ~ 

The  present  filter  is  the  top  hat  filter  and  the  closure  for  the  correlations  is  given  by 


(7) 


Lffj  +  Cali — 0 

and  the  Smagorinsky  modt'l  for  Rai  (Reynolds,  1990) 

Rai=  -  ’^I'TSai 


where  the  eddy-viscosity  is  d(>fined  by 


and 


I'T  -  (f^fi-SajSai)^ 


A  =  (Axilyi\z)3 


and 


_  I  f  dva  dvi  \ 

is  the  rate  of  strain  produced  by  the  filtered  velocity  field. 


(8) 


(9) 


(10) 


(11) 


(12) 


The  result  of  modelling  the  sub-grid-scale  motion  is  shown  in  fig.3  to  fig.9.  The 
turbulent  flow  in  the  round  jet  of  the  experiments  by  Call  and  Kennedy  (1991)  at  the 
nominal  Reynolds  number  of  Re  =  15.000  is  calculated  for  0  <  -^  <  60  without  the 
LES-model,  where  the  discretiztion  error  plays  the  role  of  the  filter  (Boris, 1990),  and  with 
the  Smagorinsky  model  discril^ed  above.  The  vorticity  magnitude  for  the  former  case  is 
shown  in  fig.3  and  in  fig.4  for  the  latter  ca,se.  It  is  evident  that  the  LES-model  using  the 
eddy  viscosity  (10)  dampens  the  smaller  scales  and  reduces  the  scale  range.  It  turns  out 
that  the  spreading  rate  in  the  LES  case  is  in  very  good  agreement  with  the  experimental 
evidence.  More  importantly,  the  effect  of  the  LES-model  on  the  particle  statistics  is  also 
beneficial  as  the  following  figures  prove.  The  dispersion  of  113/im  pentane  particles  in  the 
heated  jet  (temperature  difference  at  jet  pipe  exit:  60° K)  as  function  of  the  axial  distance 
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in  fig. 5  shows  excellnit  agreement  between  measurement  and  numerical  simulation  with 
the  LES-model.  Preliminary  results  for  the  auto-correlation  function  of  the  axial  velocity 
component  of  lighter  (oOfini  he.xadecane)  particles  in  fig.6  shows  also  close  agreement  with 
the  experiments  of  Call  and  Kennedy  (1991).  The  calculation  of  the  time  of  flight  plays  a 
noticable  role  for  dispersion  as  fig.7  illustrates.  The  experiment  is  limited  to  the  mean  time 
of  flight  wheras  the  numerical  simulationm  allows  the  computation  of  individual  times  of 
flight.  It  is  evident  from  fig.7  that  the  dispersion  calculated  with  the  experimental  (mean) 
time  of  flight  is  much  <  loser  to  the  measurements.  The  effect  of  the  LES-viscosity  on  the 
dispersion  is  shown  in  fig  8.  The  increase  in  viscosity  due  to  the  LES-model  (10)  decreases 
the  dispersion  and  improves  the  agreement  with  the  measurements.  Finally,  mean  particle 
velocity  in  axial  direction  is  fairly  independent  of  the  method  of  calculating  the  time  of 
flight  as  fig.9  shows. 
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Figure  captions 


vskip  12pt  Fig.l  Axial  mean  velocity  and  the  rms  value  of  the  axial  velocity  fluctuations 
at  several  Reynolds  numbers. 

Fig. 2  Experimental  set  uj^  for  sjnay  droplet  dispersion. 

Fig. 3  Vorticity  magnitude  in  planes  ^  =  0  and  ^  =  tt  at  dimensionless  time  t  =  898  and  a 
62  X  16  X  190  grid  for  the  flow  in  a  round  jet  (Re  —  15000)  without  the  LES-model. 

Fig. 4  Vorticity  magnittule  in  planes  ^  =  0  and  ^  =  tt  at  dimensionless  time  t  =  870  and  a 
62  X  16  X  190  grid  for  the  flow  in  a  round  jet  (Rf  =  15000)  with  the  LES-model. 

Fig.5  Dispersion  of  113///J)  vaporizing  pentane  particles  in  a  heated  jet  as  function  of  axial 
distance  obtained  with  the  LES  model  (c.,  =  0.05.  full  line)  compared  to  the  experiments 
of  Call  and  Kennedy  (1991.  broken  line). 

Fig. 6  Auto-correlation  function  for  the  axial  velocity  of  50//m  hexadecane  particles  in  an 
unheated  jet  as  function  of  the  mean  time  of  flight  (Full  line:  Numerical  simulation  with 
the  LES  model,  broken  line;  Experiment  of  Call  and  Kennedy.  1991). 

Fig. 7  Dispersion  of  113/07)  vaporizing  pentane  particles  in  a  heated  jet  as  function  of 
time.  The  numerical  results  were  obtained  w  1th  the  LES-model  and  the  time  of  flight  was 
calculated  using  the  true  Lagrangean  time  of  flight  (full  line)  and  the  mean  time  of  flight 
(broken  line)  as  in  the  experiments  of  Call  and  Kennedy  (1991.  dotted  line). 

Fig.8  Dispersion  of  113//7T?  vaporizing  pentane  particles  in  a  heated  jet  as  function  of  the 
mean  time  of  flight  obtained  without  the  LES  model  (c,  =  0.0,  full  line),  with  the  LES- 
model  (Cs  =  0.025,  broken  line)  compared  to  the  experiments  of  Call  and  Kennedy  (1991, 
dotted  line). 
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Fig. 9  Mean  axial  v<'l()(  iry  of  [)eutaue  particles  in  a  heated  je'  a.s  function  of  the 

mean  time  of  flight  obtained  with  the  LES-model  (full  line;  True  Lagrange8m  time,  broken 
line:  Mean  time  of  flight ). 


11 


siAia 


Axial  Velocity  and  Root  Mean  Square  Profiles  for 
Different  Reynolds  Number.  Nozzle  Diameter »  10mm. 
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Figure  2.  Experimental  Apparatus  for  Spray  Droplet  Dispersion 
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Appendix 


Axis  treatment  of  nonsymnietrical  flows  in  cylindrical  coordinates 

W.  Kollmaun  and  J.J.  Lienau 

MAME  Department,  University  of  California  at  Davis,  CA  95616 


1.0  Introduction 

The  turbulent  flow  in  round  jets  contains  a  rich  variety  of  structures  as  documented  in 
experiments  (Gutmark  et  ah,  1990,  Liepmann,  1991,  Yoda  et  al.,  1991,  Mungal  et  al.,  1992). 
Direct  and  LES  simulations  of  round  jets  are  thus  of  great  importance  for  the  study  of  the 
relevant  structures  in  the  development  of  the  jet  flow.  One  aspect  of  such  simulations  is  the 
treatment  of  the  flow  near  the  axis  r  =  0  in  cylindriceJ  coordinates. 

The  numerical  solution  of  the  Navier-Stokes  equations  in  cylindrical  coordinates  requires 
the  proper  treatment  of  the  discretized  equations  at  the  axis  r  =  0.  The  present  note  is 
devoted  to  this  aspect  of  the  solution  method  for  flows  that  possess  no  particuleir  symmetries. 
The  flow  domain  is  given  by  D  =  ;  0  <  r  <  Ro,0  <  9  <  25r,0  <  C  <  L)  with 

=  {{r,  d,  :  r  =  0, 0  <  ^  <  27r,0  <  (^  <  L}  as  formally  part  of  the  boundairy.  However,  it 
will  be  shown  that  there  are  no  boundary  conditions  available,  but  smoothness  and  uniqueness 
conditions  must  be  prescribed.  The  reason  for  this  is  the  fact  that  all  points  in  dT>o  are  inner 
points.  There  are  two  ways  to  deal  with  this  boundary.  The  finite  difference  grid  can  be 
staggered  such  that  r  =  0  is  not  gridline  and  no  particular  action  has  to  be  teiken,  or  the 
Navier-Stokes  equations  are  put  into  a  form  suitable  for  r  =  0  and  smoothness  and  uniqueness 
of  the  dependent  variables  are  enforced.  The  analysis  for  the  latter  case  will  be  carried  out 
and  some  results  for  the  former  case  will  be  presented. 
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2.0  Navier-Stokes  equations  in  cylindrical  coordinates 

The  method  for  the  solution  of  the  Navier-Stokes  equations  is  based  on  accurate  finite  dif¬ 
ference  schemes.  Cylindrical  coordinates  (r,  Q  appropriate  for  circular  jets  are  used.  Grid 
stretching  is  aplied  in  the  r  emd  (,*  directions  to  concentrate  the  grid  points  in  the  region 
of  interest  and  to  remove  the  outer  boundary  (Ro  denotes  the  radial  location  of  the  outer 
boundary)  as  far  as  possible  without  wasting  too  many  grid  points  (Lienau  and  Kollmann, 
1993).  However,  the  stretching  trajisformation  is  irrelevant  for  the  anedysis  of  the  variation 
of  the  dependent  variables  near  the  axis  r  =  0  and  the  standard  version  of  the  equations  will 
be  considered. 

The  Navier-Stokes  equations  for  incompressible  Newtonian  fluids  are  set  up  in  dimension¬ 
less  form  (The  jet  pipe  radius  amd  the  bulk  velocity  at  the  entrance  section  are  the  reference 
scales)  using  cylindrical  coordinates.  Mass  balance  emerges  then  in  the  form 


Id  1  dv9  dv^ 

"■^(ft’r)  ^ - ^  —  0 

r  dr  r  do  aC, 


The  radial  momentum  balance  is  given  by 


dVr  ^  V0  dVr  ^  dVr  ^  ^  1  ^  d  ^  I  d 


dur  -  .  .  „  .  . 

IT  aT  +  7  a<; 


p  dr  Re  \  dr  \  r  dr 


1  d^Vr  d^Vr  2 

■^r2  5^2  “  r2'^J 


The  azimuthal  momentum  balance  is  given  by 


dve  di'e  ve 

aT  aT 


+ 


1  d^V0  d^ve  2  dv 


-b 


r2  5^2  Q^2  ^2  QQ 


’i] 


Finally,  the  axial  momentum  balance  is  given  by 


C  uvf  ve  uvr  dvr  1  dp  I  (  1  d  /  dvr  \ 


dV(;  dv(;  Vff  dv,; 

dt  dr  ^  r  dQ 


(1] 


(2) 


'je  dve  dve  veVr  ^  dp  1  f  d  f  1  d  \ 

r  d6  5c  r  pr  d9  i?e  1 5r  \  r  dr  ) 


(3) 


1  d\^  d^v;  \ 

^r2  5^2  j 


(4) 


Inspection  of  mass  and  momentum  balances  shows  that  the  limit  r  — ►  0  produces  singular 
convective  and  viscous  terms.  Hence,  a  detailed  investigation  of  this  limit  is  called  for  to 
make  the  design  of  an  accurate  numerical  solution  procedure  possible. 
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3.0  The  limit  r  — »  0 

The  analysis  of  the  neair  eixis  variation  of  the  dependent  variables  is  based  on  the  assumption 
that  all  variables  are  at  least  three  times  continuously  differentiable  near  r  =  0.  Hence  we 
can  write 

o  q2. 

pir^e^Q)  =  p(0,^,C)  +  r^(O.^.C)  + 

and 

=  l'a(0,^,C)  +  +  ^('■^)  (6) 

valid  for  0  <  r  <  with  positive  Ca-  We  consider  first  the  scalar  variable  p  eis  r  0. 
Differentiability  at  r  =  0  implies  that  p  mtist  approach  the  same  value  as  r  — ►  0  for  any  angle 
9  and  fixed  C-  H  follows  that  the  smoothness  condition 

p(0,^,C)  =  P(0,C)  (51) 

must  hold  for  any  scalar  variable.  This  result  leads  to  the  pressure  gradient  at  r  =  0  since 
differentiation  of  (5)  with  respect  to  9  produces 

and  we  obtain 

Vp(r.«.0=  ;&(0,«C)  +5r  jg^(0,S,C)  +0(r^)  (8) 

The  smootimess  conditions  for  vectors  follow  in  aneilogous  fashion.  The  radial  velocity  com¬ 
ponent  at  r  =  0  becomes  the  angular  component  at  an  angle  9  shifted  by  and  vice  versa. 
Hence  we  have  the  smoothness  conditions 

M0,^,C,f)  =  t’<>(0,^+|,C,f)  (52) 

and 

y,(O,^,C,f)  =  «;r(O,0-^,C,O  (53) 

Furthermore,  shifting  9  hy  x  reproduces  the  radial  and  angular  components  according  to 

yr(0,^,C,f)  = + 7r,C,f)  (54) 

and 

t;9(0,^,C,f)  =  -M0,^  +  ^,C,f)  (55) 
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The  relations  (S2)  to  (S5)  are  insufficient  to  guarantuee  snxoothness  of  the  velocity  projected 
into  the  r  —  0  plane.  This  can  be  seen  as  follows.  Smoothness  implies  that  tv  and  can  be 
expanded  in  Fourier  series  at  r  =  0 


v 

iv(0,<^)  =  ^  a;exp{t(2n  -  1)^) 

fi=0 


.V 

^  a^exp(t(2n  -  1)^) 

n=0 

It  is  easy  to  see  that  (S4)  and  (So)  are  satisfied  and  (S2)  and  (S3)  lead  to 

or 

Hence,  tv  determines  vg  at  r  =  0  and  vice  versa.  The  conditions  (S2)  to  (S5)  are  thus 
recognized  as  periodicity  conditions.  Consider  now  the  vector  e(0,^.  (.f)  projected  into  the 
r  —  9  plane  as  shown  in  Fig.l.  It  is  clear  that  the  projection  is  unique.  The  projected  vector 
if  can  be  unique!  decomposed  into  the  orthogonal  components  Vr  and  vg  once  the  value  for 
the  angular  coordinate  6  has  been  chosen  (see  Fig.l). 


Fig.l.  Decomposition  of  a  vector  in  the  r  —  6  plane. 

It  follows  that 

tv(0,^,C,O  =  t’'’(C,Ocos(a'’(C,<)  -  6)  {Sr) 

and 

vg(O,0,(:,t)  =  e'’(C,Osin(a'’(C,t)  -^)  (52*) 
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hold.  It  is  clear  that  (Si*)  and  (S2*)  imply  (S2)  to  (So)  but  the  reverse  is  not  true.  Several 
useful  relations  can  be  obtained  from  (Si*)  and  (S2*)  by  differentiation.  It  follows  that 

Ov 

iv(0,^.C.O  +  ^(0.^.C,f)  =  0  (10) 

and 

y#{0.^,C,f)  -  f)  (11) 

hold. 

Finally,  the  longitudinal  velocity  component  i\-  must  approach  the  same  limit  v2Jue  as 
r  — *  0  no  matter  what  the  angular  coordinate  6  is,  hence  we  have  the  smoothness  condition 

iVd0,^.C,f)  =  ^-(0.C,f)  (S3*) 

The  following  conclusions  can  be  drawn  from  these  considerations;  The  proper  variables 
decribing  the  vector  r  at  r  =  0  are  i’P(C.f),  a^(C,f)  and  v^((,t).  The  smoothness  conditions 
(S2)  to  (So)  are  insufficient  to  determine  and  qP  uniquely. 


3.1  Mass  balance 

Application  of  the  expansion  (6)  to  (1)  leads  to 


1 

r 


M0)  +  ^(0) 


dvr  Va  dvr 


2  dr^ 


Idddr'^ 


dCdr 


(12) 


The  solution  cannot  be  smooth  as  r  — ♦  0  unless  ( 10)  holds,  which  trvkes  ceire  of  the  singularity. 
Letting  r  — »  0  we  get  then  the  proper  form  of  mass  balance  for  r  =  0 


dvr  d^t'e  dvr 


(13) 


If  differentiability  to  higher  levels  is  assumed  higher  order  versions  of  mass  balance  for  r  =  0 
can  be  obtained.  For  instaiice,  applying  (10)  and  (13)  to  (12^,  dividing  by  r  and  then  letting 
r  approach  zero  leads  to 


3  d^Vr  1  d^v$  d^vc 


(14) 


This  procedure  can  be  continued  to  arbitrary  levels  of  differentiability. 
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3.2  Radial  momentum  balance 

The  singular  part  of  the  convective  terms  in  (2)  is  given  by 


^  _  I'ef  dvr 
r[dS 

Series  expansion  of  the  expression  in  the  paranthesis  leads  to 


15) 


dVr 


dVr  , 


d^Vr 


1  a^iv 


dv 


)g  1  Vg  q 

""-o'-  ^(0)  +  O(r^ 


dd  dedr^^^'^  •l'  ar 

The  smoothness  condition  (10)  used  already  to  insure  the  nonsingulax  character  of  mass 
balance  leads  to 


dvr 

de 


-  Vg  =  r 


r 

\d0dr 


(0)-^(0)KO(r') 


proving  the  nonsingular  properties  of  convection  as  r 
given  by 


V=lfl±,  A  .  1  d^Vr  2  dvg 

’’  dr\rdr  ^  J  dd^  dd 


(16) 

0.  The  singular  viscous  terms  are 
2  dvg 


(17) 


Series  expansion  leads  to 
d^l'r  1 


d^Vr 


dvg 


VV  =  ^  +  77  ( -^(O)  +  -  2^(0)  +  ;  3^(0)  -  2 


1  /  d^Vr 


d^vg 

dddr 


(0) 


) 


The  radial  momentum  balance  should  be  nonsinguleir  as  r  =  0  is  approached  which  implies 
that  the  coefficients  of  r~^  and  r~^  must  vanish.  Fist  we  note  that  the  coefficient  of  r~^  can 
be  recast  as 


d  fdvr  \ 


Smoothness  as  r  — >  0  requires  that  (10)  and  (11)  hold,  which  implies  that  the  coefficient  of 
r“^  vanishes.  The  coefficient  of  can  be  rearranged  as 


-5^(o)-22!ii 

dS^dr  ’  dddr 


The  radial  momentum  balance  would  be  singular  if  this  coefficient  was  nonzero.  We  conclude 
that 

|^(0)-2|2(0)  =  ^(C,.)  (18) 

insures  nonsingulaxity.  The  radi2il  momentum  balance  for  r  =  0  emerges  now  in  the  form 
dVr  dVr  d^Vr  dv^  dvg  1  dp  1  f  3  d^Vr 

—  4-  J1  ,  4-  4->  tly  -  ^  f1/%  -I-  —  •  '  ' 

dt 


~  +  -  vg-^  =  - 


dr 


dddr 

d*Vr  ^  d^Vr  _  ^  &^vg  'I 


5C  dr  pdr  i2e  \  2  dr"^ 

Wd^  ~dC  ~  ~d^  } 


(19) 


strictly  valid  for  r  =  0. 
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3.3  Angular  momentum  balance 

The  singular  part  of  the  convective  terms  in  (3)  is  given  by 


_  V6  (  dve  \ 


120) 


Series  expansion  of  the  expression  in  the  paranthesis  leads  to 


Ovtheta 

dd 


dv0 


+  t’r  —  T 


I  dOdr 


(0)  +  -^(0)^  +0(r‘’) 


The  smoothness  condition  (10)  leads  to 


dvfi 

W 


-h  Vr  =  r 


I  S‘v, 
\  obOt 


(0)  +  ^(0)i>  +0(r2) 


(21) 


proving  the  nonsingular  properties  of  convection  as  r  -+  0.  The  (scalar  variable)  pressure  can 
be  expanded  and  we  get  according  to  (5) 


1  5p  _  d^p 

r  d6  dddr 

The  singular  viscous  terms  are  given  by 


(0)  +  O(r) 


d  f  I  d  \  1  d'^ve  2  dv 

7^~d9 


Series  expansion  leads  to 


V'^  = 


d'^Vg  ,  1 


dr^ 


d^Vg 


dl'r 


+  -T  -v«(0)  -f-  +  2^^(0)  +  - 


de^ 


dd 


U  ff^vg 


(22) 


r  V  d9‘^dr 


(0)  +  2 


d'^Vr 

dSdr 


(0) 


1  d^Vg  1  ^Vg  d^Vr  ^ 

The  angular  momentum  baJcince  should  be  nonsingular  as  r  =  0  is  approached  which  implies 
that  the  coefficients  of  r~^  and  r“'  must  vanish.  Fist  we  note  that  the  coefficient  of  r~^  can 
be  recast  as 

+  ^(0)  +  2^(0)  =  -MO)  +  ^(0)  +  ^  (^^(0)  +  t,.(0) 

Smoothness  as  r  ^  0  requires  that  (10)  and  (11)  hold  which  implies  that  the  coefficient  of 
r~^  vanishes.  The  coefficient  of  j — 1  can  be  rearranged  as 


d^vg 

de^dr 


d^Vr  d  (  d^Vg  dVr  \ 
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The  angular  momentum  balance  would  be  singular  if  this  coefficient  was  nonzero.  We  conclude 
that 


55^*0) +  2|^(0)- 


insures  nonsingidarity.  The  angular  momentum  balance  for  r  =  0  emerges  now  in  the  form 


dvff  dv0  d'^i'e  dv0  5iv 

dt  ^  dr  ^  ''^dedr  ac  ^  dr 


+ 


d^i's 


1  d^p 
p  d6dr 


I’ 9 


dd'^dr^'^  dC  ~dedr^j 


:24) 


strictly  valid  for  r  =  0. 


3.4  Longitudinal  momentum  balance 

The  singular  part  of  the  convective  terms  in  (4)  is  a  single  term  which  can  be  expanded  as 


(25) 


Smoothness  of  as  r  — ►  0  requires  that  i\{0.d,<^,t)  is  independent  of  the  angular  variable 
6.  hence  we  conclude  that 

ve  dv(^  d^vc^ 


r  dd 


I'e- 


holds.  The  singular  viscous  terms  are  given  by 


V  = 

r  \  dr 


:(0)  +  O(r) 

(27) 

\  1 

J  r2  de^ 

(27) 

Series  expansion  leads  to 


1 


+2  qJ  (0)  +  2  Q02Q^2  (0)  + 

The  longitudinal  momentum  baltince  should  be  nonsingular  as  r  =  0  is  approached  which 
implies  that  the  coefficients  of  and  r“*  must  vanish.  Fist  we  note  that  the  coefficient  of 
r~^  vanishes  since  r(;(0,^,C,t)  must  be  independent  of  9  according  to  (S3*).  The  coefficient 
of  must  also  vanish,  hence  we  require 


d^vc  dvr 


(28) 


to  insure  nonsingularity.  The  longitudinal  momentum  balance  for  r  =  0  emerges  now  in  the 
form 

dv,^ 


dt 


+ 


dvr  d^Vr  dvr  \  dp  1 


2dd‘^dr^  dt:^  I 


(29) 


strictly  valid  for  r  =  0. 
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3.5  Navier-Stokes  equations  at  r  =  0 

The  results  of  the  previous  section  provide  the  proper  form  of  the  Navier-Stokes  equations  for 
points  on  the  grid  surface  &Dq.  We  summarize  the  equations  first  and  then  discuss  a  possible 
strategy  for  the  solution.  Ma.ss  balance  for  r  =  0  is  according  to  section  3.1  given  by 


dvr  d^i'0  di'c 


The  momentum  balance  in  radial  direction  is  according  to  section  3.2 


dVr  dl'r 

■aT  aT 


d^Vr  dVr  dv$ 


1  r  3a^v 

p  dr  ^  Re  \  2  2 


d*Vr  d'^Vr  _  ^  d^V«  ] 

^dS^dr^  ^  dC  ~dedr^  J 

The  angular  momentum  balance  is  not  needed,  as  will  be  shown  below,  and  the  axial  mo¬ 
mentum  balance  is  according  to  section  3.4  given  by 


dv(;  ^  dv(  ^  d^t\ 

^  ^  dedr 


+  Vc 


dv(^ 


p  dC,  i2e  \  dr'^ 


1  d^l\  ) 

^2'Wd^  ] 

The  smoothness  and  uniqueness  conditions  obtained  in  section  3.0  can  be  applied  as  follows: 
The  velocity  vector  at  r  =  0  projected  into  the  r  —  9  plane  is  determined  by  and  and 
solving  the  radial  momentum  balance  at  r  =  0  for  ^  =  0  and  ^  =  f  leads  to 


Vr(0,0,C,t)  =  v^(Ct)cos(a^((,t)) 


(30) 


and 


Hence  we  get 


and 


i'r(0,^,C,0  =  y''(C-Osin(a'’(C,f)) 


0,C,<)  +  v?(0,|,C,<) 


=  arctan 


t;r(0,  f 

yr(0,0,C,O  / 


(31) 

(32) 

(33) 


for  the  projected  velocity  at  r  =  0.  Once  and  have  been  determined  the  radial  and 
angular  velocity  components  at  r  =  0  follow  from  (Si*)  and  (S2*). 

The  strategy  for  the  solution  includes  the  following  steps  to  deal  with  the  axis  r  =  0: 
Solve  mass  balance  ( 13),  the  axial  momentum  balance  (29)  and  the  radial  momentum  balance 
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(19)  for  ^  =  0  and  d  =  j.  Then  we  obtain  the  projected  velocity  vector  from  (32)  and  (33). 
The  radial  and  circumferential  velocity  components  follow  then  from  (Si*)  and  (S2*)  as 


Vr(0.0,i;j) 


i’;^(0,0,C,O  +  i’r(0- C- 0 ‘-‘0!>(arctan 


rv(0,0,C,O  / 


(34) 


and  _ 

re(0,^,C,t)  =  yi>2(0,0,C,<)  +  t’r  ( 1 1 C  ^  0  sin(  arct  an  ^  (35) 

This  completes  the  treatment  of  the  axis  if  r  =  0  is  gridline. 


4.0  Application 

The  finite-difference  method  described  by  Lienau  and  Kollmann  (1993)  was  applied  to  the 
prediction  of  the  flow  in  a  round  jet  for  a  nominal  Reynolds  number  of  Re  =  15000.  The 
development  of  the  axial  and  circumferential  instabilities  leads  quickly  to  solutions  without 
symmetries,  hence  requirung  the  proper  treatment  of  the  dependent  variables  near  the  axis 
r  =  0.  Two  treatments  of  the  axis  are  considered: 

(1)  The  grid  contains  r  =  0  as  gridline  and  the  smoothness  conditions  (Si)  to  (S5)  are 
enforced  but  not  the  uniqueness  conditions  (Si*)  to  (S3*). 

(2)  The  grid  is  staggered  and  r  =  0  is  not  gridline. 

The  results  for  case  ( 1 )  in  Fig.2  for  velocity  and  in  Fig. 3  for  vorticity  show  that  spurious  higher 
harmonics  in  circumferential  direction  appear  and  are  insufficiently  damped.  The  solution 
remains  bounded  but  the  near  eixis  variation  of  velocity  and  vorticity  are  unsatisfactory.  The 
results  for  the  case  (2)  in  Fig. 4  for  velocity  and  in  Fig. 5  for  vorticity  on  the  other  hand  are 
quite  satisfactory. 


5.0  Conclusions. 

A  finite  difference  method  for  the  Navier-Stokes  equations  in  cylindriceil  coordinates  was 
developed  that  is  sixth  order  accurate  in  space  and  second  order  accurate  in  time.  The 
treatment  of  the  dependent  variables  and  the  structure  of  the  Navier-Stokes  near  r  =  0  was 
established  using  series  expansio;\s  and  the  conditions  that  variables  and  equations  must 
remain  nonsingular  and  unique  for  -►  0.  The  application  of  two  different  treatments  of  the 
axis  r  =  0  lead  to  the  conclusion  that  staggered  grids  are  superior  over  non-staggered  grids, 
where  the  enforcement  of  smoothness  is  not  sufficient  to  produce  satisfactory  results.  It  is 
clear  that  uniqueness  conditions  should  be  imposed  in  this  case,  but  this  has  not  yet  been 
tested. 
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Fig. 2  Velocity  in  planes  0  =  0  and  0  =  rr  at  dimensionless  time  t  =  4C  for  the  How  in  a  round 
jet  (Re  =  15000)  for  case  (1). 
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Fig.5  Voitifity  in  planes  6*  =  0  anti  6  -  n  at  diniensionless  time  t  =  57  f(»i  tin*  How  in  a  ronnil 
jet  (Re  =  15000)  for  etise  (2). 


